DNA methylome and transcriptome profiling reveal key electrophysiology and immune dysregulation in hypertrophic cardiomyopathy

ABSTRACT Hypertrophic cardiomyopathy (HCM) is the most common inherited heart disease. However, a detailed DNA methylation (DNAme) landscape has not yet been elucidated. Our study combined DNAme and transcriptome profiles for HCM myocardium and identify aberrant DNAme associated with altered myocardial function in HCM. The transcription of methylation-related genes did not significantly differ between HCM and normal myocardium. Nevertheless, the former had an altered DNAme profile compared with the latter. The hypermethylated and hypomethylated sites in HCM tissues had chromosomal distributions and functional enrichment of correlated genes differing from those of their normal tissue counterparts. The GO analysis of network underlying the genes correlated with DNAme alteration and differentially expressed genes (DEGs) shows functional clusters centred on immune cell function and muscle system processes. In KEGG analysis, only the calcium signalling pathway was enriched either by the genes correlated with changes in DNAme or DEGs. The protein-protein interactions (PPI) underlying the genes altered at both the DNAme and transcriptional highlighted two important functional clusters. One of these was related to the immune response and had the estrogen receptor-encoding ESR1 gene as its node. The other cluster comprised cardiac electrophysiology-related genes. Intelliectin-1 (ITLN1), a component of the innate immune system, was transcriptionally downregulated in HCM and had a hypermethylated site within 1500 bp upstream of the ITLN1 transcription start site. Estimates of immune infiltration demonstrated a relative decline in immune cell population diversity in HCM. A combination of DNAme and transcriptome profiles may help identify and develop new therapeutic targets for HCM.


Background
Hypertrophic cardiomyopathy (HCM) is a common hereditary cardiovascular disease with a global incidence of~0.2% [1,2]. It is a major cause of sudden death and its main clinical manifestation is unexplained left ventricular hypertrophy (≥15 mm). Physiological characteristics of HCM include hypercontraction of the myo fasciculus, impaired diastole, increased sensitivity to calcium ions, and perturbed metabolic signalling [3,4]. Changes in immune status may also occur during the tissue pathogenesis of HCM [5]. Immunocytes account for 5.3% of the total number of ventricular cells [6]. They maintain cardiac homoeostasis and normal physiological function but also mediate adverse post-injury inflammatory responses, myocardial remodelling, and cardiac electrophysiological abnormalities [7][8][9][10]. Cardiac hypertrophy may cause immune dysregulation that negatively impacts cardiac function [7]. HCM patients are at a relatively high risk of sudden death, and heart transplantation is the only efficacious treatment. However, transplantation is often impeded by long wait times because of the shortage of histocompatible donors and the risk of immune rejection. Hence, innovative therapeutic modalities are urgently required to delay the progression of HCM and mitigate the risk of sudden cardiac death associated with it.
DNA methylation (DNAme) is a type of epigenetic modification that affects the entire human genome. It enables the same DNA sequence to be transcribed differently in various cell types. DNAme regulates tissue and organ growth and development. DNAme abnormalities are implicated in the pathogenesis and histopathological changes related to various diseases [11]. CpG hypermethylation in promoter regions downregulates gene transcription. By contrast, CpG hypermethylation in the gene body indicates active gene transcription [12][13][14]. DNA methyltransferases (DNMTs) and ten-eleven translocation cytosine dioxygenases (TETs) regulate DNAme and demethylation, respectively [15,16].
Aberrant DNAme is involved in the growth and pathological adaptation of cardiomyocytes [17,18]. DNMT3A knockdown by CRISPR/Cas9 creates aberrant DNAme in cardiomyocytes which, in turn, modulates the expression of genes associated with altered myocardial contractility, damages the mitochondria, and causes defects in lipid and glucose metabolism [19]. Changes in the DNAme alter the genes encoding potassium channels which play critical roles in cardiac conduction and arrhythmia [20]. Increasing evidence supports the involvement of DNAme in cardiomyopathy, heart failure and transplantation, immune effects and distant effects (e.g., myocardial remodelling and activation of fibrotic processes) [21][22][23][24][25][26][27][28][29][30]. In heart failure patients, for instance, satellite (SAT) element hypermethylation is associated with significant upregulation of SAT transcription and may play a role in maintaining chromosomal integrity and myocardial stress response [21,22]. DNAme predicts the risk of complications after transplantation, and FOXP3 promoter hypermethylation predicts immune rejection of heart transplantation [23,24]. In late-stage HCM, some patients may have heart failure. Thus, studies investigating heart failure may detect DNAme alterations occurring in HCM. Nevertheless, dynamic changes in DNAme exist across distinct stages in disease progression [28,31,32], and the profile for HCM patients prior to heart failure remains unexplored. Furthermore, aberrations of DNAme may regulate innate and adaptive immunity [33][34][35]. However, little is known about the mechanisms of immune infiltration and aberrant DNAme in HCM myocardium.
We performed a microarray-based genomewide DNAme analysis paired with transcriptome data to compare HCM against normal left ventricle myocardium. The combination of DNAme screening and transcriptome signals enabled us to investigate HCM mechanisms. It also permitted us to evaluate the effects of the previously hypothesized roles of immune cell function, cardiomyocyte development, and cardiac electrophysiology on HCM. The potential of the immune-related genes identified herein for the treatment of HCM could be explored in subsequent functional studies.

Myocardial specimens
All subjects provided written informed consent prior to participation. All procedures were approved by the Ethics Committee of Beijing Anzhen Hospital (No. 2021155×). Thirty-two myocardial tissue samples were collected from 24 patients with HCM and eight patients with normal myocardia at Beijing Anzhen Hospital. The HCM samples were obtained during septal myectomy and transported to the laboratory in ice-cold cardioplegia solution until cryopreservation (<4 h). HCM diagnosis was defined according to the 2020 AHA/ACC Guideline for the Diagnosis and Treatment of Patients with Hypertrophic Cardiomyopathy [36]. All patients were diagnosed by two experienced clinicians. HCM was confirmed by 2D echocardiography showing unexplained left ventricular hypertrophy with diastolic interventricular septal thickness≥15 mm and septal wall: posterior wall thickness ratio≥1.3 in the absence of any other cardiac or systemic disease causing a similar magnitude of hypertrophy. Normal, healthy myocardial tissues were collected from donor hearts of patients who voluntarily donated their bodies for research. For all donors, clinical examination and medical history displayed no indications of cardiac history nor structural heart disease. All hearts were arrested in situ using ice-cold cardioplegia solution and transported to the lab on wet ice (always<4 hours). The left ventricles (LVs) of hearts were dissected and flash frozen in liquid nitrogen for this study. Tissue DNA was isolated with a QIAmp DNA Mini Kit (Qiagen, Hilden, Germany) and quantified by Qubit fluorometer (Thermo Fisher Scientific, Waltham, MA, USA). Total DNA>0.5 ug was selected and all samples were used for DNAme data collection. Total RNA was isolated from 300 mg tissue with TRIzolTM reagent (Invitrogen, USA) according to the manufacturer's instruction. RNA concentration was measured using Qubit® RNA Assay Kit in Qubit® 2.0 Fluorometer (Life Technologies, CA, USA). RNA integrity was checked using the RNA Nano 6000 Assay Kit on the Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA). The selection criteria were total RNA>0.4 μg, and RNA integrity number (RIN) > 5. Of these, 23 (16 from HCM patients and 7 from normal controls) were qualified for further sequencing.

RNA-seq protocol
Poly(A) RNA was purified from 0.4 μg total RNA using Dynabeads Oligo(dT)25-61005 (Thermo Fisher, USA). The poly(A) RNA was fragmented into small pieces using Magnesium RNA Fragmentation Module (NEB, USA). The cleaved RNA fragments were reverse-transcribed into cDNA by SuperScript II Reverse Transcriptase (Invitrogen), which was subsequently used to synthesize U-labelled second-stranded DNAs with E. coli DNA polymerase I (NEB), RNase H (NEB) and dUTP Solution (Thermo Fisher). Size selection was performed with AMPure XP beads (Beckman Coulter, USA). After the heat-labile UDG enzyme (NEB) treatment of the U-labelled second-stranded DNAs, the ligated products were amplified with PCR. The average insert size for the final cDNA library was 300 ± 100 bp. The 2 × 150 bp pair-end sequencing for mRNA were performed on Illumina Hiseq 2500 (Genechem Technology Co. Ltd., Shanghai, China).

DNA methylation
DNA was bisulphite-converted with the EZ DNAm kit (Zymo Research) according to the manufacturer's instructions and then hybridized to the Infinium MethylationEPIC BeadChip (850 K, Illumina). These microarrays were scanned using the Illumina HiScan SQ scanner by Emei Tongde Ltd.

Histology
HCM and control LV tissue samples were successively perfused by saline and 4% paraformaldehyde. LV tissues were fixed for 24 h (4% paraformaldehyde), transferred to ethanol (70%) for subsequent dehydration and paraffin embedding. LV samples were continuously cut into sections, and stained with Masson's trichrome and Hematoxylin-eosin (HE) staining by following the protocols.

DNAme data QC and differentially methylated probe (DMP) analysis
The bisulphite treatment in the EZ DNA Methylation-Gold Kit was used to convert the DNA samples according to the manufacturer's instructions (Zymo Research, Irvine, CA, USA). DNAme was detected with an Infinium Human Methylation EPIC BeadChip (850k) Microarray (Illumina, San Diego, CA, USA). Raw fluorescence data for the DNAme were stored in IDAT files. The 'ChAMP' package (https://bioconductor.org/ packages/release/bioc/html/ChAMP.html) in R (R Core Team, Vienna, Austria) was used for QC and data analysis. Raw data were loaded with the 'champ.load' function and transformed into βvalues. The latter represent the ratios of the fluorescence intensities of methylated probes to those of unmethylated probes and range from 0 to 1. The β-value increases with the degree of methylation. This function was also used for initial probe screening and QC. Probes meeting any of the following criteria were removed: (1) p > 0.01; (2) < 3 beads in≥5% of all samples; (3) non-CpG; (4) multihits; (5) underlying SNPs; and (6) location on X or Y chromosome. The filtered probe set was used in the subsequent analyses. Beta mixture quantile expansion (BMIQ) normalization and singular value decomposition (SVD) were used to eliminate batch effects caused by differences in experimental time points. A principal component analysis (PCA) was performed using the 'prcomp' function in the 'stats' package (https://statistics globe.com/stats-r-package) of R. PCA reduced the dimensionality of the high-dimensional methylation site matrix/gene expression matrix into a form with only a few principal components. This small number of principal components can then effectively represent the variation at all methylation sites/gene expression, and therefore describes the differences between samples.
The 'Bumphunter' method (https://bioconductor. org/packages/release/bioc/html/bumphunter.html) was used to screen DMPs. The △β-value is the arithmetic difference between the mean β-value of the HCM group and that of the normal control for a single probe. Correction of P values (adjusted P) for multiple comparisons was performed by using the Benjamini-Hochberg procedure. The DMP thresholds were adjusted P < 0.05 and log|FC| > 0.1. The latter term is equivalent to the effect of |△β-value| > 0.1. HCM hyper DMPs had logFC>0.1 whereas HCM hypo DMPs had logFC<0.1. DMP gene annotations were performed using the default annotation of the 'ChAMP' package in R. The reference genome was hg19 (https://www.ncbi.nlm.nih. gov/assembly/GCF_000001405.13/). An in-house Python script was used to evaluate probe overrepresentation in each chromosome. To determine DMP overrepresentation in a chromosome in all DMPs (DMPSet), it was necessary to determine the a priori distribution of the number of probes (N) from a chromosome in a random probe set (ProbeSet R ) equal in size (M) to the DMPSet. To this end, M probes were sampled 10,000 times with replacements from all probes in the array, and the cumulative probability of observing N > M within each chromosome in a ProbeSet R was determined. For HCM hyper DMPs and HCM hypo DMPs, 2,202 and 1,864 probes, respectively, were randomly sampled 10,000 times from 732,724 total probes in the array. The significance of the overrepresented chromosome was the sum of the probabilities (P-values) of finding N > M probes from this chromosome in the DMPSet. P value<0.05 indicated that the number of probes was significantly increased on the chromosome.
Those with adjusted P < 0.05 and log 2 |FC| > 1.0 were considered differentially expressed genes (DEGs). PCA was performed using the built-in function in DESeq2.

Gene set enrichment analysis (GSEA) and gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses
GSEA and GO and KEGG enrichment analyses were performed using the 'clusterProfiler' package (https://bioconductor.org/packages/release/bioc/ html/clusterProfiler.html) in R [40]. All DMPs were annotated with hg19 to obtain the DMP genes. The GO (adjusted P < 0.05) and KEGG(P value<0.05) enrichment analyses were performed on all DMP genes and DEGs. The significantly enriched GO items were grouped in either dataset into the functional network within ClueGO v. 2.5.9 (https://apps. cytoscape.org/apps/cluego) in Cytoscape v. 3.9.1 (https://github.com/cytoscape/cytoscape/releases). The default parameters were edge-weighted, forcedirected, and BioLayout (http://biolayout.org/down load.html) for CluePedia (https://apps.cytoscape. org/apps/cluepedia). Cytoscape was used to construct a network from the overlapping pathways in the results of the GO enrichment analyses of the DMP genes and the DEGs. The κ (kappa) score threshold was 0.3 and the sharing group percentage was 30.0%. GSEA was applied only to DEGs with adjusted P < 0.05.

Protein-protein interaction (PPI) network analysis and pathogenic gene screening for overlapping genes of DMP genes and DEGs
A PPI network of overlapping genes was constructed based on the STRING (https://string-db.org/), an online protein network database [41]. Genes with moderate confidence interaction scores (>0.4) were selected and isolated node genes were removed. The remaining genes were clustered into subnetworks by the Markov clustering (MCL) method. The inflation parameter was 1.5. All other parameters used in this analysis were set by default. The immune-related gene dataset that was confirmed in earlier studies to be implicated in the natural immune process was downloaded from the IMMPORT database [42] (https:// www.immport.org/shared/). The GeneCards (https:// www.genecards.org/) database was used to screen and obtain pathogenic HCM genes. The keyword was 'hypertrophic cardiomyopathy' [43].

xCell immune infiltration analysis
The xCell analysis uses a deconvolution algorithm and integrates the advantages of GSEA. The xCell calculated single-sample gene set enrichment analysis (ssGSEA) scores for 489 gene signatures. Gain compensation was then corrected and new scores were used to assess the relative abundance of 64 cell types in each tissue sample. Cell types included multiple adaptive and innate immune, haematopoietic progenitor, epithelial, and extracellular stromal cells as well as 48 TME-associated cells. Gene expression values (fragments per kilobase per million mapped reads or FPKM) from the bulk RNA-seq data were uploaded to http://xCell.ucsf.edu/ [44] to obtain 64 cell scores per sample. Immune cells with P value<0.05 (t-test) were selected for display.

Data availability
The raw DNAme data (No. PRJCA009134; https:// ngdc.cncb.ac.cn/omix/preview/ocAxwKJt) and the raw RNA-seq data (No. PRJCA009145; https:// ngdc.cncb.ac.cn/gsa-human/s/mRwWy2L3) reported herein were deposited to the Genome Sequence Archive [45] of the China National Center for Bioinformation/Beijing Institute of Genomics of the Chinese Academy of Sciences and is publicly accessible.

Clinical characteristics of study patients and donors
Thirty-two cardiac ventricular septal tissue samples were collected from 24 HCM patients and eight donors (Methods; Table 1). DNAme profiles were determined for all samples, and paired RNA-seq analyses were performed on 23 samples (Figure 1a). The HCM group consisted of 14 males and 10 females with average age 56 ± 13 y and maximum left ventricular wall thickness (LVST) = 19.3 ± 4.4 cm. There was no significant difference between the sexes in terms of their maximum LVST (P = 0.6854; t-test). HE staining results indicated myocyte hypertrophy, myocyte disarray, nuclear hyperchromatism, and vacuolar degeneration in HCM patients compared with normal controls in 40× and 200× views. Interstitial fibrosis and replacement fibrosis were detected in Masson-stained samples of HCM patients (Figure 1b). The HCM patients harboured different variants in the exons of their pathogenic HCM genes MYBPC3, MYH7, MYL2, MYL3, TNNI3, TNNT2, TPM1, and ACTC1 [2]. No major pathogenic variants were responsible for patient stratification (Methods; Supplementary Table S1). We did not find significant differences in either methylome or transcriptome in terms of gender (male, or female), BMI (>24 or<24) and diabetes (with or without).

DMPs: methylation level and chromosome distribution
After quality control, 732,724 CpG probes were obtained and used to cluster 32 samples into two groups consistent with disease grouping (PCA; Methods; Figure 1c). The HCM myocardium presented with both hypomethylation and hypermethylation changes and exhibited elevated double β-value distribution peaks (Methods; Figure 1d). There were 4,066 DMPs (adjusted P < 0.05; log|FC| > 0.1; Methods; Supplementary    Figure 1e). These were −0.29 to −0.1 for the HCM hypo DMPs and 0.1 to 0.34 for the HCM hyper DMPs. The DMPs were unevenly distributed across the chromosomes (Genomic annotation; Methods). The HCM hyper DMPs were enriched on chromosomes 4,5, and 12 while the HCM hypo DMPs were located on chromosomes 2, 3, 9, 10, 13, 18, and 21 (P value<0.05, Figure 1h). In terms of genomic distribution (Figure 1f,g), 48% and 42% of the HCM hypo DMPs and HCM hyper DMPs, respectively, were located within the gene body while 84% and 66% of HCM hypo DMPs and HCM hyper DMPs, respectively, were located within the opensea. Overall, the methylation levels of the 4,066 DMPs significantly differed and could effectively discriminate the HCM myocardium from the normal tissues (Figure 1i).
A GO enrichment of the DEGs and DMP genes revealed pathways with extensive DNAme and gene expression alterations that were functionally grouped into various networks (Methods; Figure 2i). Two major immune subnetworks included most of the GO items. One was characterized by 'regulation of leukocyte migration' while the other was characterized by 'muscle system process.' Calcium signalling was the only significantly enriched KEGG pathway for both the DMP genes (P value = 2.18 × 10 −4 ) and the DEGs (P value = 6.04 × 10 −3 ) (Figure 2j).
The transcriptional activity of the DNMTs and TETs did not explain the observed changes in HCM DNAme. There were no significant differences in transcription between HCM myocardium and normal tissues for any DNMT family member genes (DNMT1, DNMT3A, and DNMT3B; P value>0.05; t-test; Figure 2e) or TET-related genes (TET1, TET2, and TET3; P value>0.05; t-test; Figure 2g). A public HCM dataset validated our observations (P value>0.05; t-test; Figure 2f,h). Of the methyl-CpG-binding domain (MBD) proteins or epigenome readers [46], EGR1 was downregulated and BAZ2B was upregulated, but there were no DNAme alterations (Supplementary Table  S4). Changes in DNAme were detected in MBD4, ZBTB38, KLF4, WT1, and UHRF2 but without any obvious alteration of transcription (Supplementary  Table S2). Compared with previous observations in human cardiomyocytes (hCMs) [47], DMP genes influenced all three constructed function groups: muscle contraction, cardiac transcriptional regulation, and heart development (Supplementary Table S5).

The nodes of the PPI network of DEGs harbouring DMPs underscored the roles of immune regulation, cardiac development, and electrophysiology in HCM pathogenesis
To identify DEGs caused by DNAme alterations, we extracted 85 genes common to both the DEGs and the DMP genes (Supplementary Table 6; Figure 3a, left). Twenty-four percent (20/85) of them were HCM genes in GeneCards (Supplementary Table 7). TTN-AS1 was the most relevant gene (score = 169) to HCM. Hence, these disease-causing genes may also have DNAme alterations and contribute to HCM occurrence. Nine of the 85 common genes, namely, ESR1, GDF6, IL20RA, PMP2, SCG2, STC2, TGFA, TGFB2, and VEGFC, are associated with cytokines, cytokine receptors, the TGFb family, interleukin receptors, and antimicrobials (IMMPORT database; Supplementary Table 8). PPI networks were resolved for 79 DEGs with aberrant DNAme (Methods; Figure 3a, right). Genes associated with disease causality, namely, ESR1, RUNX2, and CACNA1A, were crucial nodes with edges>3. ESR1 and RUNX2 belonged to the larger of the two clusters associated with tissue and organ growth and development. CACNA1A was the central node of the cluster involved in cardiac electrophysiology.

Dysregulated immune responses were observed in cardiac hypertrophy and remodelling
The heart consists of heterogeneous cell subtypes intensively involved in cardiac hypertrophy and remodelling [6,75]. We estimated the immunocyte composition (macrophages, T cells, and dendritic cells (DCs)) of the HCM myocardium using the novel gene signature-based xCell method (xCell; Methods). The ImmuneScore was significantly lower for the ventricle tissues of the HCM group than it was for those of the normal group. Macrophages comprise a wide range of functionally heterogeneous phenotypes. Understanding the roles of various macrophage phenotypes during cardiac remodelling could help develop a promising therapeutic strategy [7]. In the present study, we found that the total macrophages (P value<0.01; t-test) and M1 macrophages (P value<0.01; t-test) were significantly decreased in HCM patients (Figure 3c). Macrophages respond to different environmental signals and activate various polarization programs [76,77]. Similarly, monocytes (P value<0.05; t-test) and induced dendritic cells (iDCs; P value<0.05; t-test) were significantly decreased in HCM patients. Relative to the public dataset, our results showed that the active dendritic cells (aDCs) (P value<0.01; t-test) and M2 macrophages (P value<0.05; t-test) were significantly decreased in HCM myocardial tissue (Methods; Figure 3d). The preceding results reflect end-stage HCM which might differ from early disease onset.
The stromal cells were significantly increased in the public datasets (GSE130036) as well as our own compared with the normal group (Figure 3c,d). The foregoing findings indicate that nonmyocytes in general and immune cells, in particular, were dysregulated in HCM patients. The immune response and regeneration capacity may be significantly impaired in HCM patients as their immune, stromal, and stem cells are all dysregulated. Nevertheless, the mechanisms of these dysregulations remain to be elucidated and merit further investigation.

Discussion
DNAme alterations have been detected in various diseases and are promising as diagnostic and prognostic biomarkers [78]. As an inherited disease, abnormal myocardial function due to pathogenic genomic alterations is anticipated to be an important aetiological factor of HCM, and we believe that altered DNA methylation may be a secondary alteration in myocardial function compensation. Recent progress has shown methylation changes may result from altered transcription activity in response to a stimulus in cardiomyocytes [79]. The authors found that the acute hypoxia stress response continuously activated specific gene expression patterns and resulted in DNAme changes in regulatory regions of the corresponding genes. Therefore, transcriptional and DNAme changes help elucidate these molecular mechanisms of HCM pathogenesis. However, it requires considerable effort to determine DNAme profiles for human cardiac tissues. Animal models of heart failure and tissue samples from heart failure patients showed DNAme alterations in targeted genomic regions [27][28][29], here we further provided DNAme profiles from HCM patients in surgery prior to heart failure to explore changes in a relatively earlier stage in disease progression. Unlike tumour cells, HCM cardiomyocytes retain their basic morphology and functions. As transcriptional changes in methylation-related enzymes do not occur in HCM myocardium, there may be fewer changes in HCM DNAme than in the tumour genome. Abnormal gene expression in tumours may explain the extensive epigenetic reprogramming characteristic of genome-wide DNAme alterations. Future research should endeavour to investigate whether changes in DNMT and TET protein abundance occur in human HCM myocardium and identify their contributions to alterations in HCM DNAme. Besides, methylation quantitative trait loci (meQTL)/ Expression quantitative trait loci (eQTL) analysis was used to correlate genetic variant loci with DNAme levels at specific CpG loci/expression of specific genes to explain genetic variation in the disease [80,81].
With accumulation of this type of data for HCM, meQTL/eQTL analysis will provide a better explanation of genetic variation in familial HCM.
Myocardial hypercontraction, dysregulation of calcium homoeostasis, and metabolic signalling disorders are pathological manifestations of HCM [4,[82][83][84]. To the best of our knowledge, the present study provided DNAme profiles for HCM tissues with paired transcriptome datasets and determine the alterations in DNAme underlying changes in gene expression. DNAme profiles facilitate the exploration of the connections underlying myocardial remodelling and cardiac electrophysiological abnormalities during pathogenesis. They may also disclose candidate targets for the development of novel therapeutic strategies for HCM.
We observed decreases in innate and adaptive immune activity as well as mononuclear macrophage system-related components in the HCM myocardium. These changes are closely related to cardiac electrophysiological abnormalities and myocardial remodelling. Macrophage depletion leads to abnormalities in cardiac electrical signalling and especially atrioventricular block and predisposes the heart to progressive cardiomyopathy, reduced cardiac output, diastolic dysfunction, and impaired haemodynamics [10,85]. The most significant change in HCM myocardium was a reduction in iDC content. These antigenpresenting cells have a strong migratory capacity, and various types of iDCs reside in the  [86]. In dilated cardiomyopathy, apoptosis and insufficient angiogenesis decrease DC diversity and increase the number of mature DCs [87]. This mechanism may explain the reductions in iDCs observed in HCM myocardium. In comparison, DEGs in hCMs were enriched only in pathways associated with cell lineage and tissue differences rather than immune-related ones [47]. Immune cell population changes in HCM are likely to result from myocardial remodelling or myocardial injury, but further efforts are required to dissect their impacts on cardiomyocytes. Future integration of high-resolution single-cell methylation sequencing and single-cell transcriptome sequencing may reveal methylation alterations in immune cell and cardiomyocyte interactions. We found that ESR1 encodes a cytokine receptor also as a transcription factor that has motifs for estrogen binding, DNA binding, and transcription activation. It had the most edges in the HCM PPI network. The two HCM hyper DMPs within its gene body may explain its upregulated transcription level. Overexpression of the ESR1 protein product ERα might be a protective factor in disease progression, prevent cardiac hypertrophy [61,62,88], reduce arrhythmias [69], improve vascular endothelial function [89], and modulate innate immune signalling in DCs and macrophages [49]. There were no differences between the sexes in terms of their myocardial ESR1 expression levels [90][91][92]. Nevertheless, males have relatively higher incidence and significantly lower average age at HCM diagnosis than females [93][94][95][96]. The hormone receptor ESR1 might help maintain basic physiological functions in lesions and slow disease progression. Future investigations into the therapeutic mechanism of ESR1 might be warranted.
ITLN1 underwent the most significant transcriptional changes and corresponding DNAme alterations in all HCM samples. In HCM myocardium, ITLN1 has a hypermethylated site within 1500 bp upstream of the transcriptional start site, which may lead to its transcriptional downregulation. ITLN1 is a component of the innate immune system. It is an anti-inflammatory adipocytokine present in epicardial adipose tissue and it has cardioprotective efficacy [50]. A recent study demonstrated altered transcriptional activity in failing ventricles [97]. Our HCM patients displayed significantly downregulated ITLN1 and, therefore, limited macrophage differentiation into the anti-inflammatory M2 phenotype [51]. We also observed decreased macrophage composition in the transcriptome dissection for our HCM samples. Hence, ITLN1 downregulation is a key factor in this phenomenon. Although the present work characterized HCM at the DNA and RNA levels, these bulk tissue data can not identify changes specific to a cell population, and future efforts are required to determine the DNAme alterations in cardiocytes and validate the composition changes of immune cells in HCM myocardium.
In the present study, we profiled the changes in DNAme in HCM myocardium and evaluated their corresponding impact by jointly analysing DNAme and the transcriptome. We built upon traditional pathological alterations in HCM and identified possible links among immune dysregulation, cardiac electrophysiological abnormalities, and myocardial remodelling. Immunerelated genes that also regulate DNAme might help identify and develop novel therapeutic strategies for HCM. cncb.ac.cn/omix/preview/ocAxwKJt while the raw RNA-seq data are hyperlinked to the dataset at https://ngdc.cncb.ac.cn/ gsa-human/s/mRwWy2L3.

Consent for publication
Not applicable.

Disclosure statement
No potential conflict of interest was reported by the authors.

Ethics approval and consent to participate
The experimental protocol was established according to the ethical guidelines of the Helsinki Declaration and was approved by the Ethics Committee of the Beijing Anzhen Hospital, Capital Medical University (No. 2021155×). Written informed consent was obtained from individual participants or their legal guardians.